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Summary 

Skelton’s  original  Fortran  program  for  predicting  time-harmonic  sound 
radiation  from  an  infinite  cylinder  comprising  elastic,  viscous  and 
acoustic  layers  has  point  forces  and  monopoles  as  excitations.  It  has  been 
extended  to  cover  the  case  of  excitation  by  a  plane  wave  at  arbitrary 
incidence,  and  the  output  may  be  used  for  transient  radiation/scattering 
studies.  Numerical  examples  are  given,  but  a  computer  listing  is  not 
included. 
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FORTRAN  PROGRAM  FOR  SOUND  SCATTERING 
BY  LAYERED  ISOTROPIC  CYLINDER 


By  J  H  James 
1 .  INTRODUCTION 

Theoretical  aspects  of  acoustic  scattering  from  infinite  cylindrical 
geometry,  insonified  at  an  oblique  angle,  have  been  given  by  Flax  et  al  [1] 
for  the  case  of  a  solid  elastic  cylinder  and  James  [2]  for  the  case  of  a 
fluid  filled  cylindrical  shell.  Skelton  [3]  has  given  the  theory  necessary 
to  calculate  sound  radiation  from  a  layered  cylinder  excited  by  time- 
harmonic  mechanical  forces  or  acoustic  monopoles,  the  layers  being  composed 
of  acoustic  fluids,  viscous  fluids  and  elastic  solids.  Skelton's  analysis 
and  computer  program  have  since  been  modified  to  cover  the  important  case 
of  sound  scattering  due  to  plane  wave  excitation  at  oblique  incidence,  and 
this  substantially  forms  the  content  of  this  report.  Theoretical  aspects 
are  not  discussed  as  they  can  be  found  in  the  references  cited. 

2.1.  General 

Figure  1  shows  a  cylindrically  layered  system  comprising  isotropic  elastic 
solids,  viscous  fluid  layers  and  acoustic  fluid  layers;  governed  by  the 
exact  linear  equations  of  elastodynamics ,  viscodynamics  and  acoustics.  The 
time-harmonic  excitations  are:  mechanical  force  at  elastic  layer  interface; 
acoustic  monopole  in  exterior  fluid;  acoustic  monopole  in  interior  fluid; 
acoustic  monopole  in  acoustic  layer;  and  acoustic  plane  wave  incident  at  an 
oblique  angle. 

The  scope  of  the  Fortran  program,  which  calculates  far-field  acoustics,  is 
evident  from  its  specification  which  can  be  found  in  the  Appendix. 

2.2.  Time-Harmonic  MF  Scattering 

In  order  to  demonstrate  the  capability  of  the  program,  consider  a  sandwich 
construction  comprising  facings  of  concentric  steel  cylinders,  and  a  much 
softer  filler  material.  The  SI  constraints  of  the  3  layers,  LamS  material 
properties  A  and  p  being  in  GN/m2,  are  taken  as: 


Material 

A 

P 

V 

a 

b 

Top  Facing 

104.4 

75.6 

7700 

0.01 

0.99 

1.00 

Soft  Filler 

0.04 

0.01 

900 

0.01 

0.91 

0.99 

Bottom  Facing 

104.4 

75.6 

7700 

0.01 

0.90 

0.91 

The  exterior  water  has  sound  velocity  '500  m/s  and  density  1000  kg/m3,  and 
there  is  a  vacuum  inside  the  inner  steel  cylinder,  p  is  material  density, 

>7  is  the  hysteretic  loss  factor,  a  is  inner  radius  and  b  is  outer  radius. 
Figure  2  shows  monostatic  target  strengths  (see  Appendix  for  definition)  of 


a  cylinder  comprising  the  outer  steel  facing  alone.  The  frequency  range 
extends  to  the  medium  acoustic  ka  value  of  8.4.  The  features  have  been 
analyzed  by  Veksler  [4]  and  others  in  terms  of  the  first  symmetric  Lamb 
mode  (S0)  and  the  first  antisymmetric  (Aq)  Lamb  mode,  with  the  former  mode 
dominating.  In  a  flat  plate  these  modes  correspond  to  the  familiar 
longitudinal  and  flexural  waves.  Figure  3  shows  corresponding  target 
strengths  of  the  composite  cylinder.  In  the  frequency  range  to  800  Hz, 
interactions  between  the  facings  via  the  stiffness  of  the  filler  are 
especially  evident  at  normal  incidence.  Above  this  frequency,  the  outer 
facing  becomes  more  and  more  decoupled  from  the  inner  facing. 

2.3.  Time  Harmonic  HF  Scattering 

Figure  4a  shows  monostatic  target  strengths,  at  normal  incidence,  of  a 
steel  cylinder  of  outer  radius  1  m  and  thickness  2  cm.  The  frequency  range 
extends  to  the  high  acoustic  ka  value  of  104.7.  The  dips  below  10  kHz  and 
rapid  oscillations  above  this  frequency  are  predominantly  due  to  the 
Lamb  S0  and  Aq  modes  respectively,  the  dips  being  explained  in  terms  of 
resonances  and  the  rapid  oscillations  as  due  to  coincidence  lobe  radiation. 
Figure  4b  shows  the  corresponding  plot  when  the  thickness  of  the  cylinder 
is  increased  to  10  cm.  The  influence  of  higher  order  Lamb  and  surface  wave 
modes  must  be  present  to  a  certain  extent,  but  physical  explanations  in  the 
published  literature  are  unclear.  Both  these  plots  resemble  the 
corresponding  plots  obtained  for  spherical  geometries  by  White  [5],  but  the 
oscillations  herein  are  much  less  marked  due  to  weaker  resonant  behaviour 
by  cylinders. 

2.4.  Transient  S.ca.t  feting 

Time-harmonic  excitation  has  been  assumed  throughout,  but  transient 
excitations  are  readily  available  using  a  recently  developed  software  tool 
[6]  for  transient  response  calculations  from  time-harmonic  spectra.  For 
example,  Figure  5a  shows  the  (total)  far-field  sound  level  spectrum,  at 
9=90  and  4=0,  of  the  composite  cylinder  when  a  time-harmonic  monopole  is 
located  close  to  its  surface.  Figure  5b  shows  the  transient  far-field 
sound  pressure  when  the  excitation  is  a  400  Hz  sine  wave  which  is  on  for  a 
single  cycle.  The  early  response  (note  the  time  scale  is  relative)  is 
dominated  by  direct  radiation  from  the  monopole  and  to  a  lesser  extent  by 
its  specular  reflection  from  the  cylinder  surface.  The  extended  ringing  in 
the  tail  is  of  small  amplitude. 

3.  CQB£m.SIQM§ 

The  numerical  results  demonstrate  that  the  Fortran  program,  arising  from 
Skelton's  theoretical  work  extended  to  cover  the  case  of  excitation  by  a 
plane  wave  at  oblique  incidence,  is  a  useful  research  tool  for  analyzing 
the  acoustics  of  cylinders  comprising  isotropic  layers.  Time-harmonic 
responses  are  the  standard  output  of  the  program,  but  transient  responses 
are  readily  available  using  auxiliary  software. 

4.  ACajQWlEPSBffiMI 

The  not  inconsiderable  effort  of  E  A  Skelton  when  doing  the  original 
theoretical  work  and  computer  programming  is  much  appreciated. 
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APPENDIX  A 


Program  Specifications  of  CYL001.FQR 

Sound  radiation  and  scattering  by  infinite  layered  cylinder,  see  Figure  1 
for  geometry,  comprising  isotropic  layers  of  elastic  solid,  viscous  fluid 
and  acoustic  fluid. 

Item  01 : 

N  The  number  of  layers. 

Item  02:  for  each  layer,  in  outer  to  inner  order,  give  the  following. 

LAM  real  Lame  first  constant. 

MU  real  Lame  second  constant. 

ELAM  hysteretic  loss-factor  for  first  Lame  constant. 

EMU  hysteretic  loss-factor  for  second  Lame  constant. 

RHO  density. 

B  outer  radius  of  layer. 

A  inner  radius  of  layer. 


according  to  the  table: 


LAM 

b.d 

ELAM 

EMU 

RHO 

0.0 

0.0 

Vc 

c 

P 

for  acoustic  layer 

0.0 

M 

0.0 

c 

P 

for  viscous  layer. 

A 

<n 

n* 

P 

for  elastic  layer. 

in  which  c  is  sound  velocity;  p  is  density;  A  and  p  are  Lame  constants  for 
the  elastic  layer,  and  p  is  the  dynamic  coefficient  of  viscosity  for  the 
viscous  fj.uid  layer,  r)^,  r\u  and  rjc  are  hysteretic  loss- factors ,  viz. 
c=c(i-irjc),  A=A(1-ifjj),  p=p(  1  -irju)  .  Note  a  fluid  loss  factor  is 
equivalent  to  a  structural  loss  factor  of  2r;c.  If  an  elastic  solid  has 
Young's  modulus  E  and  Poisson's  ratio  v,  then  the  equivalent  Lame  constants 
are  A=Eu/(  1 +u)  ( 1 -2u)  and  p  =  E/2(1+u). 

Item  03: 


RHOEXT  density  exterior  acoustic  fluid. 

CEXT  sound  velocity. 

Item  04: 

RHOINT  density  interior  acoustic  fluid.  Set  to  0.0  if  vacuum. 

CINT  sound  velocity.  Set  to  0.0  if  vacuum. 


Item  05: 


ITYPE  =0  for  point  force  excitation. 
=1  for  exterior  point  source. 
=2  for  interior  source. 
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=3  for  source  inside  layer. 

=4  for  plane  wave  excitation. 

Item  06:  only  if  ITYPE=0 

NF  number  of  point  forces,  then  for  each  force: 

FCOM  amplitude  of  force;  complex  number,  eg  (1.0, 0.0). 

IFACE  interface  number  of  point  of  application,  the  outer  surface  being 
interface  1,  etc. 

IRZ  =  1/2/3  for  force  acting  in  radial/tangential/axial  direction. 

ZF  axial  coordinate  of  force. 

Item  06:  only  if  I TYPE  =  1,2,3 

NS  number  of  point  sources,  then  for  each  source: 

CAMP  amplitude  of  source;  complex  number,  eg  (1. 0,0.0). 

ILA  layer  number  of  source  location,  the  outer  layer  being  layer  1, 

etc.  Set  to  0  for  ITYPE  =  1,2. 

RS  radial  coordinate  of  source. 

ZS  axial  coordinate  of  source. 

Item  07 : 

THETA  elevation  angle  in  degrees  of  observation  point. 

PHI  circumferential  angle  in  degrees  of  observation  point. 

Item  08: 

NFREQ  number  of  frequencies. 

FMIN  start  frequency. 

FMAX  final  frequency. 

Item  09: 

NBES0  terms  in  Fourier  summation  at  first  frequency. 

NBES1  terms  in  Fourier  summation  at  last  frequency. 

1.  For  ITYPE  =  0,1, 2, 3  the  program  output  is  a  table  of  far-field  sound 
pressure  level  in  dB  ref.  1  micropascal  at  1  m,  assuming  (if  logic  permits 
the  user  to  do  so)  that  the  excitation  amplitudes  are  normalized  to  rms,  ie 
an  excitation  of  (1.0, 0.0)  is  assumed  to  be  an  rms  amplitude  of  unity.  The 
first  column  of  output  is  the  frequency;  the  second  column  of  output  is  the 
far-field  sound  level;  the  third  column  of  output  is  the  same  as  the  second 
column  for  ITYPE=0,2,3  and  for  ITYPE=1  it  is  the  far-field  pressure 
assuming  a  rigid  immovable  exterior  boundary.  Excitations  for 

ITYPE=0, 1 , 2 , 3  are  applied  at  the  circumferential  angle  d= 0.  The  angles  9 
and  ^  of  the  far-field  observation  point  are  selected  in  Item  07. 

2.  For  ITYPE=4  the  second  column  of  the  program  output  is  the  target 
strength  of  the  elastic  system  and  the  third  column  is  the  target  strength 
of  the  system  assumed  rigid  and  immovable.  The  incident  plane  wave  for 
ITYPE=4  comes  from  the  direction  ^=180.  The  angle  of  measurement  i  should 
be  chosen  as  180  for  monostatic  target  strength.  Several  hours  of 
contemplation  will  reveal  that  for  this  substantially  two-dimensional 
scattering  problem,  the  azimuthal  angle  of  incidence  0,  is  identical  to  the 
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azimuthal  angle  $  of  the  observation. 


The  target  strength  is  defined  as 


T  =  20.  log10[(r/Pi)1/2|psf(r)tf.  ,*,*..=180)1] 


where  Px  is  the  amplitude  of  the  incident  wave  and  ps£  is  the  pressure 
scattered  to  the  far-field. 

3.  Input  is  from  the  file  FOROOS.DAT  and  output  is  to  the  file 
FOR006.DAT.  Results  are  also  put  on  disk  file  in  a  format  suitable  for 
plotting  by  a  system  dependent  program  GRAPH1 . 

4.  Maximum  values,  which  can  be  changed  by  altering  the  Fortran  PARAMETER 
statement  are:  NLAY=1 0 ,  NFREQ=501 ,  NF=10,  NS=10,  NBES1=100 

5.  Field  quantities  of  interest  are  expanded  as  Fourier  series  in  the 
circumferential  coordinate  Initially,  run  the  job  with  five  frequencies 
and  select  and  change  NBESO  and  NBES1  until  the  answers  appear  to  have 
converged  to  a  sufficient  accuracy.  Start  with  NBES0=5  and  NBES1=5+ka, 
say.  Note  Bessel  functions  in  the  Fourier  series  may  cause  overflow  if 
NBESO  or  N8ES1  is  too  large.  Ideally  the  program  should  be  run  in  a 
Fortran  environment  in  which  64  bit  arithmetic  with  a  large  exponent  is 
used . 


6.  Fatal  error  numbers  in  output: 

10  =  error  when  calculating  system  inverse  in  LAYER. 

11  =  error  when  calculating  layer  matrix  inverse  in  SCALC. 

7.  An  auxiliary  Fortran  program  CYL002  is  available  for  calculating  real 
wavenumber  versus  frequency  plots  for  selected  circumferential  harmonics. 
These  dispersion  plots  help  in  the  physical  interpretation  of 
radiation/scattering  spectra. 
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TARGET  STRENGTH 


FIG. 2  Monostatic  target  strength  of 
thin  steel  cylinder,  (a)  8=90  (b)  0=50 
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FIG. 3  Monostatic  target  strength  of 
composite  cylinder,  (a)  0=90  (b)  9=50 
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FIG. 5  Far-field  radiation  from 
close  to  composite  cylinder,  (a) 
harmonic  spectrum  (b)  transient 
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